Scale-dependent halo bias from scale-dependent growth 
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We derive a general expression for the large-scale halo bias, in theories with a scale-dependent 
linear growth, using the excursion set formalism. Such theories include modified gravity models, 
and models in which the dark energy clustering is non-negligible. A scale dependence is imprinted 
in both the formation and evolved biases by the scale-dependent growth. Mergers are accounted for 
in our derivation, which thus extends earlier work which focused on passive evolution. There is a 
simple analytic form for the bias for those theories in which the nonlinear collapse of perturbations 
is approximately the same as in general relativity. As an illustration, we apply our results to a 
simple Yukawa modification of gravity, and use SDSS measurements of the clustering of luminous 
red galaxies to constrain the theory's parameters. 
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I. INTRODUCTION 

The cosmological constant (A) + cold dark matter 
(CDM) -I- general relativity (GR) model has been very 
successful in accounting for current cosmological data. 
It is incumbent upon us to put this standard model to 
further tests. One of its characteristic predictions is a 
scale-independent sub-Hubble linear growth, which can 
be violated if dark energy clusters^ , or if gravity is mod- 
ified [2HS]. Both types of models introduce new scales 
into the growth of structure: the Jeans scale in the case 
of clustered dark energy, and the GR-to-non-GR tran- 
sition scale in the case of modified gravity. The most 
direct way to test for this effect is to measure the matter 
power spectrum at different redshifts, and reconstruct the 
growth factor as a function of scale. Here, we focus on a 
corollary of a scale-dependent growth, a scale-dependent 
halo bias, which in principle allows us to discern scale 
dependence even with measurements of the large-scale 
structure at a single redshift. 

The fact that a scale dependence in growth implies 
scale dependence in halo bias (on linear scales) was 
pointed out by [TD] (henceforth HP), generalizing earlier 
work by [nilI2]: 



D{ko,Zio 



(1) 



D{ko, Zobs) 

where 61(^0? ^obs; ^form) signifies the linear bias on scale 
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^ Dark energy other than the cosmological constant generally clus- 
ters, though the degree of sub-Hubble clustering is often negli- 
gible. Exceptions include models where the dark energy sound 
speed is substantially less than unity, e.g. [l]. 



ko observed at redshift Zobs, for haloes that form at red- 
shift Zform- The symbol D denotes the linear growth fac- 
tor at the relevant scale and redshift. This expression, 
which assumes passive evolution, i.e. halo number con- 
servation after formation, tells us that the observed bias 
would inherit a scale dependence from the growth, even if 
the formation bias is scale-independent. In this paper, we 
wish to relax these two assumptions: scale-independent 
formation bias, and no mergers. 

The paper is organized as follows. The extended 
Press-Schechter, or excursion set, formalism [I3l - fl7] is 
described and generalized to allow for a scale-dependent 
growth in Sec. |lTj This is used to compute the halo mass 
function and the halo bias. In Sec. |III[ we describe an 
illustrative example, a modified gravity model of the 
Yukawa type, and present a calculation of the linear 
growth factor and the excursion barrier (collapse thresh- 
old). We present the halo bias for this example, and 
compare it with observations. We conclude in Sec. |IV| 
In this paper, for the purpose of illustration, we have cho- 
sen to focus on one particular model of modified gravity. 
Our results in Sjll] for the scale-dependent halo bias are, 
on the other hand, fairly general. A recipe for using 
these results in more general settings is summarized in 



Sec. IV Readers who are interested primarily in applica- 
tions, and not on the derivation, can skip directly to Sec. 
|IV[ In Appendix [X] we discuss some details concerning 
the use of the excursion set method for theories with a 
scale-dependent growth factor. In Appendix [B] we give 
the derivation of the Yukawa model from a scalar-tensor 
theory, and we describe our spherical collapse model in 
Appendix [Cj 

Before we proceed, let us briefly discuss the connection 
with some of the literature on the subject. The halo mass 
function for a Yukawa theory like the one we study was 
computed by [18]. Our paper follows their formalism, and 
it is in a sense a straightforward extension to compute 
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the conditional mass function and halo bias. Halo bias 
in f{R) gravity and the DGP models has been measured 
from numerical simulations in P^^E5] . They find fair 
agreement between simulations and the bias derived from 
a modified Sheth and Tornicn [24^ mass function whose 
parameters reflect the altered spherical collapse. 

There is a large literature on testing GR using the 
growth of large-scale structure, e.g. [25]. Most studies 
allow for a scale-dependent growth factor, but ignore its 
effect on the galaxy bias. Our expression for the bias 
should be useful for incorporating the latter effect into 
such studies. 

There is also a substantial literature on a large-scale 
scale-dependent bias from primordial non-Gaussianity, 
e.g. [IB]. This other source of scale dependence should 
be distinguishable from that from growth; we will discuss 
this in Sec. IVl 



II. EXCURSION SET THEORY OF HALOES 
WITH A SCALE-DEPENDENT GROWTH 
FACTOR 

We describe in detail here the excursion set formalism. 
Much of the discussion replicates the standard treatment, 
but with special care taken to allow for a scale-dependent 
growth factor. 



A. Random walks 

The premise of the excursion set theory is that a halo 
will form from a region of a certain size when the over- 
density of matter S, smoothed over that region, is greater 
than a critical value Jcrit [IS]- In FIG.[lj the mass over- 
density, at a point, smoothed on a comoving scale i?, is 
plotted against the variance S on that scale around that 
point, where 

P(fc) is the matter power spectrum, and W{kR) the 
Fourier-space window function. As the smoothing scale 
is decreased, the smoothed density field undergoes a 
random walk. The halo mass M is given by M — 
(47r/3)p(z — 0)R'^, where p is the mean matter density 
today. Hence, M, S and R can all be thought of as 
equivalent variables. 

Following common practice, we smooth the density 
field with a real-space top-hat filter (even though, strictly 
speaking, a true random walk requires a top-hat filter 
in Fourier space |15|V In an important sense all of our 
analysis is 'native' to real space; for all quantities not 
explicitly marked as being in Fourier space by a ~, we de- 
fine the wavenumber associated with a scale R as simply 
k = 1/R. The Fourier transform of a real-space top hat 




large R: large mass S(R.) —Cp" ^i^^" ^- small mass 
in region ^ ' in region 



FIG. 1: Random walks: notation for the excursion set model. 
A point is contained in a halo corresponding to the first up- 
crossing of the barrier, here (5crit, by the random walk. 

is 

3 

WikR) — {sin(a;) — x cos(x)} , x = kR . 

Following |18j , this random walk exercise is carried out 
at some initial redshift Zi, chosen to be early enough 
that the linear growth is scale-independent (i.e. neither 
dark energy nor modified gravity is important by this 
redshift^). Throughout this paper we set Zi = 100, al- 
though its precise value is not important as long as it is 
sufficiently early. We add the subscript i to the critical 
threshold for halo formation (5crit,i to remind ourselves 
that it is the comparison between the smoothed dm at Zi 
against this threshold that defines haloes. The value for 
<5crit,i in general depends on both the halo mass of inter- 
est (i.e. the barrier is not necessarily flat, unlike in the 
standard model) , and the actual physical redshift of halo 
formation Zform- 

B. Unconditional first-crossing distribution &: the 
mass function 

The central quantity in the excursion set theory is the 
first-crossing distribution, /. Define /(S*!, (5i|0, 0) dS*! as 
the probability that a random walk starting at the ori- 
gin will first cross the barrier, of height Si = (5crit,i, 
at 5* = S'l ± dS'i/2. This is equivalent to saying that 
/(S'i,5i|0,0)|dS'i/dMi|dMi is the probability that our 



^ In the standard model of GR + A or (nearly) homogeneous 
dark energy, the linear growth is scale-independent, the redshift 
at which the density field is defined and the random walk is 
carried out is immaterial. This is discussed in Appendix [A] 
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point (around which this smoothing or random walk exer- 
cise is performed) belongs to a halo of mass Afi ±dMi/2. 
We can think of {Si,6i) as the target point we are aiming 
for on our random walk. 

To find /, it is sometimes convenient to consider first 
the cumulative crossing distribution, F, which is the 
probability that the walk first crosses the barrier at 
S < Si: 



F{Si,5i\0,0)= ['\ 
Jo 



f{S,S\0,0)dS, 



where S is the critical threshold (5crit,i, generally a func- 
tion of S; in the case where it is S'- independent, i.e. flat 
barrier, F (and therefore / = dF/dSi) can be found us- 
ing a clever trick due to Chandrasekhar j27j . or by mod- 
elling the random walk as a diffusive process (see P^^llSj ). 
For a more general barrier, an alternative method for 
finding / will be discussed in |IIF| We assume that the 
probability distribution of the density field as a function 
of the smoothing scale R is Gaussian: 



PG{6,SiR))dd = 



exp [-5'^/2S{R)] d5 . 



Using the above, we see that the mass function of 
haloes (the physical number density of haloes as a func- 
tion of mass, smoothed over the whole universe) which 
are identified at Zform and observed at Zohs is given by 



n{Mi)dMi 



p{Zohs) 

Ml 



/(^i, -5110,0) 



ASi 



dMi 



dMi , (3) 



where p(2:obs) is the mean matter density at z^hs- The 
mass function has an implicit dependence on Zform 
through the dependence of the threshold/barrier 5i on 
the formation redshift. 



Conditional first-crossing distribution & halo 
bias 



1. Lagrangian space 

We can now move the starting point of our random 
walk away from the origin, to some other point on the 
S-d plane. Call this point (5*0, (5o)- This is equivalent 
to calculating all our quantities within a region of the 
size on which the mass variance is Sq, and which has a 
mass overdensity 6q. The differential probability that a 
random walk starting from this point will first cross the 
barrier Si is f{Si^ 5i\Sq, Sq). In the flat-barrier case this 
only depends on 5*1 — Sq and Si ~ Sq. 

Using this conditional first-crossing distribution we can 
construct an expression for the overdensity of haloes 
within a region of {Sq,So) [28l |29]- Just as we derived 
the unconditional mass function above, we can derive a 
conditional mass function. We can see from our defini- 
tions that the mean number Af of Mi haloes in a region 



of total mass Mq is 



U{Mi,Si\Mq,Sq) 



Mq 

Ml 



f{Sl,Sl\SQ,SQ) 



dSi 



dMi 



(4) 



and so their number density is M/Vq if the region has 
volume Vq{Mq). 

We define the overdensity of haloes 'I' in a region '0' 
as the fractional overabundance of these haloes in this 
region compared to their number density smoothed over 
the whole universe, giving 



^^(1|0) = 



N{Mi,Si\Mo,Sq) 
n{Mi,Si)VQ 



- I 



(5) 



This is the halo overdensity in Lagrangian co-ordinates, 
in which Mq « p{zi)VQ, defining Vq. Indeed, all quantities 
in this Lagrangian space discussion, including 5q,Vq,5i 
and so on, are defined at the initial redshift Zi. 



2. Eulenan space 

In Eulerian co-ordinates we need to replace Vq with the 
evolved Eulerian volume. In the above, the Lagrangian 
space quantities Vq and Sq are calculated at z^. We now 
have three redshifts: the initial redshift Zi (which we al- 
ways set to too), the identification redshift of the haloes 
•Zform (in a sense, their formation redshift), and the red- 
shift at which the haloes are observed Zobs- 

Label as S the physical matter overdensity of the region 
of interest, at the observation redshift Zobs- It will be 
related to the Lagrangian space Sq by 



D{Sq, Zobs) r 

D[So,Zi) 



(6) 



because Sq is defined at z^. The use of the linear growth 
factor D is justified for a small 5*0, or correspondingly 
a large region. By the time of observation, the vol- 
ume of the observed patch will have contracted from its 
Lagrangian-space value, because conservation of mass re- 
quires that the volume decreases as S increases. Label as 
V the volume of our region at Zi , taking into account the 
contraction it will experience by the time it reaches Zobs- 
Since S{zi) « 0, we find 



Mq^p{z,)Vo = p(z,){1 + S)V 



V = 



Vq 



Substituting Vq V into Eq. ([5|, and using Eqs. (H 
and H, gives the Eulerian halo overdensity Sj^: 



5f: 



(l + '5) 



/(5i,,5i|^o,^o(^)) 
/(5i,5i|0,0) 



(7) 



It is clear that St (dropping the 'E' from now on) will 
in general not be the same as the matter overdensity in 
the same region, S. A positive S will usually help haloes 
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reach the collapse threshold and increase their overden- 
sity, with the opposite effect for negative i5. There are 
also halo exclusion effects; for example, the halo density 
will be suppressed unless their mass is much smaller than 
the mass available in the region. From the random walk 
perspective, the 'time' available for barrier crossing is 
Si — Sq, in which time the walk must travel a 'distance' 
6i~do- You can suppress halo formation either by reduc- 
ing the time available (Mhaio = Mi approaching Mq) or 
by increasing the required distance (increasing the height 
of the barrier). 

Since we know Sh is a function of 6, we can expand it 
as a Taylor series [30j . 



OO ^ 



(8) 



fc=0 



The bk co-efRcients are the bias parameters. The first, 69, 
vanishes in the limit of Sq <^ Si. Here, we are primarily 
concerned with 61, the linear bias, in the (5 — ?> limit (be- 
cause we are interested in large regions, on which scales 
the mass overdensity is expected to be small) . The linear 
bias in the limit (5 — > is just dSh/dS evaluated at (5 = 0, 
and so, using d6o/dd = D{Sq, Zi) / D{Sq, Zohs), we find 



6i((5^0) 



/(Si,'5i|0,0) 



f{Si,5i\So,0) 



D{So,z,) df{Si,Si\So,So) 
D{So,Zohs) d6o 



(9) 



Recall that 5*0, Rq, and kg are all equivalent variables; 
So and Rq are related by Eq. ([2]), and fco = l/^o- The 
D{So, Zi)/ D{So, Zohs) factor in this equation will make 
the large-scale bias scale-dependent in many modified 
gravity theories, or theories in which dark energy clus- 
ters in a non- negligible manner. 



D. Flat barrier approximation 



As we will see in pil[ a flat barrier may be a good ap- 
proximation for some theories. This is not surprising even 
for non-standard theories: haloes of a sufhciently small 
mass cross the barrier in a regime where essentially the 
standard story applies, i.e. the crossing occurs on small 
scales where gravity is Newtonian and dark energy typi- 
cally is quite smooth. A useful analytic form for the halo 
bias can be written down in this case. The conditional 
first crossing distribution is jl5j 



f{Si,Si\So,do) = 



1 



Si 



/2n {Si - 5o)3/2 



exp 



[Si - So)' 



2{Si - So) 



and so, using Eq. ^, we find the Eulerian-space bias to 
be 



bi = 



Si 

Si ~ Sq 

X exp 



3/2 



D{So,Zi) ( Si 



1 



D{So, Zobs) \Si - So Si 



SlSo 



2Si{Si — So) 



(10) 



We are interested in scales on which So ^ 5*1, i.e. 
the halo mass of interest is much smaller than the mass 
encompassed in the region over which we are calculating 
the clustering. Therefore, 



D{So,Zi) f 61 



^1(5*07 ^obs; ^foim) — l+T^TT^, ' \ [ ^ 1 

D(So, Zobs) V Si 



1 

Ti 



(11) 



The bias bi depends on Zform because i5i is the initial 
overdensity threshold for collapse at Zform- 

As a check, it can be seen this reduces to the 
Mo & White result in the standard model, by 
recognizing that both Si and 5*1 in our descrip- 
tion are defined at the initial redshift Zi — one 
can use part of the factor D{So, Zi)/ D{So, Zoi^s) = 
[D{zi)/D{zform)][D{zioTm)/D{zohs)] to rcscale Si^/Si and 
Si down to Zform (recalling that D is scale-independent 
in the standard model), and obtain bi — 1 + 
[D{ziorm)/D{zohs)]{'^'^ - l)/^c, wherc Sc is the linearly 
extrapolated overdensity of collapse at Zform'^, and v = 
Sc/ \/ Si, form with Sijorm being the the variance at Zform 
for our halo mass of interest Mi. The bias in its most 
familiar form obtains when one sets Zform = -Zobs- 

In non-standard models, E 
bias 61 is scale (So) dependent in general 



(11) implies the linear 
The scale 

dependence comes entirely from the ratio of growth factor 
-D(S'o, 2;i)/D(S'o, Zobs)- The result resembles Eq. ([T]) [ID] 
To facilitate comparison, it is helpful to rewrite Eq. ( |11[ ) 
as 



6l(<So, Zobs; Zform) = 1 + (&1 (5*0, Zform) " 1) 



D{So, Zform) 



D{So,Zohs) ' 

(12) 

where &i(5'o, Zform) is the formation bias, given by 
h \ ^ I D{So,z,) fSi^ \ 1 

bi{So,ziorm) = 1 + vh^"M:r (^3) 

L'(5o, Zform) \Si J di 

This is a trivial rewriting, but it makes clear that 
the formation bias 61(50, 2form) is in general scale- 
dependent, contrary to what was assumed in HP. 
It also helps differentiate between two different ef- 
fects: one is the scale dependence at formation (from 



^ For instance, 5c ~ 1.686 for a universe with a critical matter 
density and no cosmological constant or dark energy. More gen- 
erally, for a flat universe with matter and A, 5c ~ 1.686 X (1 -|- 

0.123 logiof^/), with Uf = !^m(l + Zform)^/Pm(l + ^form)^-|-f^A] 
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D{So, Zi)/D{So, Zform) in Eq. [Is]) , and the other is the 
scale dependence from passive evohition thereafter (from 
D{So,Zforni)/D{So,Zohs) in Eq. [12]). 

Note that when Zform 7^ ^obs, the haloes are identified 

at 2;form ^ ^obs ; 

and maintain their identities until Zobs ■ 
Mergers in the form of accretion onto these haloes are 
allowed between the two epochs, so long as the identifi- 
cation is unaltered, i.e. the haloes, once identified, are 
conserved. If this is not true — if the haloes continue 
to merge with each other all the way to Zobs — then one 
should use Xform = -Sobs- In this case, there is no distinc- 
tion between the bias observed at Zobs and the formation 
bias. 

It is also instructive to further rewrite the formation 
bias using the following relations: 



— <5l,forr 
Si — S'l form 



DjSuz,) 

D{Si, Zform) ' 

DHSi,z,) 



(14a) 
(14b) 



where we have introduced ^i^form and S'l^form for the lin- 
ear overdensity threshold and the variance at Zform (recall 
that our random walks are performed at the initial red- 
shift Zi, hence 61 and 5*1 are defined then). The formation 
bias can thus be written as 



bi{So, Zfo 



X , £'(5i,Zforni) [ ^ 



form 



1 



) \s 



l.form 



Si 



foiT 



(15) 

where we have assumed that at a sufficiently early Zi, 
the growth factor is scale-independent, i.e. D{SQ,Zi) — 
D{Si,Zi). This is a useful expression because <5i^form 
is typically a constant, with fairly weak dependence on 
Zform and cosmology (see footnote [3| . This appears to 
be a good approximation even for the modified gravity 
model we will study below, where for sufficiently small 
haloes 5i,form approaches its ACDM value, which is 1.671 

for Zform = 0. 



E. The modified linear growth factor 

So far, we have not been very explicit about where the 
scale-dependent growth factor actually comes from. On 
the sub-Hubble scales of interest, mass and momentum 
conservation imply the following equation for the linear 
matter overdensity S: 



d^J 
da2 



1 dH\ dS _ 1 



— = 



H da J da a^H^ 



(16) 



* These expressions work even for 



Zfo 



< 2obs- This is relevant if 



one identifies haloes at some low redshift Zformi but is interested 
in the clustering of the center of mass of their constituents at 
some earUer redshift z^hs- 



where the Laplacian is in comoving coordinates. What 
modified gravity, or clustering in dark energy, does is 
to modify the relation between S and the gravitational 
potential (j), from the standard Poisson equation. The re- 
sult can often be modelled as replacing hy O ■ 5, 
where O is some linear operator which can be scale- 
dependent, giving rise to a scale-dependent growth for 
d. Solving the resulting equation in Fourier space would 
give 6{k,z) = [D{k, z)/ D{k, Zi)]S{k, Zi), where k is the 
wavenumber of interest. We use the symbol D to dis- 
tinguish this growth factor from the growth factor D 
we use in the excursion set calculation. The difference 
arises from the fact that the excursion set calculation 
cares about the growth of the variance S, smoothed on 
scale R, or equivalently, M. We define D in such a way 
to give the correct evolution of the variance S: 



D^{R,z) 



S{R,z) 



S{R,Zi) 



(17a) 



(2vr)3 



P(fc,zO <^ W{kR) 



D{k,z) 
b{k,Zi) 



(2^)^ 



P{k,z,)W^{kR). (17b) 



In the excursion set computation, the overdensity S is 
evolved using the same growth factor D. In practice, the 
integral over k is expected to be dominated by fc ^ 1 /R, 
and therefore D{R,z) is roughly D{k = 1/R,z), though 
we do not use this approximation. 



F. Calculating the first-crossing distributions 

As mentioned above, the flat barrier approximation 
turns out to be a fairly good one for the example we 
study. However, in general, the barrier or threshold 
for collapse is mass/scale-dependent. The corresponding 
first-crossing distribution can be found as follows. 

First of all, it is possible to calculate the / distribution 
analytically for barriers which are linear functions of S: 
if Si {S)=uj~l3S, 



f{Si,Si\So,5o) 



AS 



exp ■ 



{AS - pAsy 

2AS 



where AS = SiiSo) ~ So and AS = Si - So [11]. For a 
more general barrier, there is no exact analytic solution, 
but reasonable approximations have been worked out by 
Sheth & Tormen, and Lam & Sheth [351 [33] . In this pa- 
per, we adopt the algorithm of Zhang & Hui [34] which 
gives an exact, albeit numerical, solution for the uncon- 
ditional /-function. It is straightforward to extend their 
method to the conditional case. It can be shown that the 
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function is the solution of the integral equation 

f{Si,6i{Si)\So,Sa) = gi[Si,6i{Si); So,6o] 

rSi 

dS[ g2[Si,Si{Si); S[,di{S[) 



where 



X f{S[,Si{S[)\So,So), 



01 — OQ QOi 



g2[Si,Si{S,);Sl6,{S[)] 



diSi Si — tS^ 

X Pg{5i{Si)--5i{S[),Si-S[). 



This is a Volterra equation of the second kind. By 
treating / as a vector /, we can see that this equation is 
of the form / = g+Mf where g is another vector and M 
is a matrix, implying that / can be found with a matrix 
inversion: f — [I — M)~^g. M is a triangular matrix, 
allowing I~M to be inverted efhciently by iteration [55] . 



III. AN ILLUSTRATION WITH A YUKAWA 
MODEL 

By the Yukawa model, we mean using the following 
effective Poisson equation to compute the growth of per- 
turbations and the formation of haloes (for instance, it 
goes into the right hand side of Eq. [I6] ) : 



1 + a 



1 



(18) 



which is written in Fourier space. Here, Gjv is the effec- 
tive Newton's constant at high k. The k"^ /{k"^ + a^/A^) 
term corresponds to a Yukawa modification, where 1/A 
can be thought of as a mass: it determines on what scale 
modifications to GR become important. The parameter 
a quantifies the size of the modification, i.e. the ratio 
of effective Newton's constant on small scales to on large 
scales is 1/(1 -I- a). A positive a corresponds to weaken- 
ing gravity on large scales (large compared to A), while 
a negative a (but > — 1) corresponds to the opposite. 
This effective Poisson equation can be solved by 



(r) 



G 



N 



1 



'r' f^ {l + ae^^{-\r-r'\/\)} 
\r — r'\ 

(19) 

where r is the physical proper coordinate, related to the 
comoving coordinate x hy r = ax. 



The model encapsulated by Eq. ( 18 ) is phenomenolog- 



ical in nature. For a positive a, it can be derived from a 
scalar-tensor theory with a potential for the scalar field, 
described in Appendix [B| 



One implication of Eq. ( 19 1 is that Birkoff 's theorem, 
or Gauss's law, no longer holds: a test particle in a 



spherically-symmetric matter distribution feels a force 
from the matter both inside and outside the surface of 
constant radius defined by its position. 

We investigate three methods, one exact and two ap- 
proximate, for solving for the nonlinear evolution of 
spherically-symmetric perturbations, the details of which 
are deferred to Appendix [C] This nonlinear problem has 
to be solved to determine the appropriate collapse thresh- 
old (5crit,i used in the excursion set calculation. The 
violation of Birkoff's theorem implies an initially top- 
hat perturbation does not remain so. However, we find 
that modelling the spherical perturbation as homoge- 
neous (i.e. a top hat) at every time step is a good 
approximation for the parameters and halo masses of 
interest. Other infra-red modifications of gravity, such 
as DGP, appear to share the same feature — Schaefer & 
Koyama have argued [36] that the effect of the violation 
of Birkoff's theorem in DGP gravity should be small. 

We find that a Yukawa-type modification of gravity, 
which is consistent with observations of galaxy clusters 
and hence has A > 5 Mpc, will only affect the collapse 
dynamics for haloes with mass M > 10^'^h~^MQ; for 
smaller haloes the evolution is basically GR-like. This is 
consistent with previous studies of spherical collapse of 
large haloes in Yukawa-modified gravity [18] . 

In the rest of this section, we present results for the bar- 
rier (5crit,z from the spherical collapse computation, and 
for the halo bias 61 from the excursion set method. The 
background cosmology is taken to be that of ACDM (see 
Appendix |B] for discussions on consistency with modify- 
ing gravity for the perturbations). For the power spec- 
trum, we use the analytic fit from [37] . with parameters 
fixed by the WMAP 5-year analysis: flcDM = 0.206, 

^^baryon — 0.043, il\ = 1 — (f^CDM + ^^baryon) , "-s — 0.961, 

h — 0.724 f35]. We normalize the power spectrum at a 
high redshift (zi = 100), i.e. the primordial normaliza- 
tion is the same for both the standard and non-standard 
models. The amplitude is chosen such that, if GR were 
to hold, (78 would equal 0.787 today. 



A. Barriers 

The barrier used in the random walk calculation is the 
initial overdensity required for a perturbation to collapse 
to zero radius by an identification redshift Zform- In gen- 
eral relativity the required initial overdensity is indepen- 
dent of scale, i.e. the barrier is flat. Curved barriers have 
been used to model the effect of ellipsoidal, rather than 
spherical, collapse; e.g. see [55] . 

For the Yukawa model, the barrier can be generated 
using the spherical collapse calculation outlined in Ap- 
pendixjC] in what follows we use the homogeneous sphere 
approximation throughout. For each scale ki (= l/i?i, 
or a mass of Mi), the initial overdensity is found which 
collapses at Zform- This defines (5crit,i, which is equated 
with Si in the excursion set theory, as a function of fci or 
Ml. FIG. [2] shows barriers for the Yukawa model with 
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FIG. 2: Barriers in Yukawa models parametrized by a and A: 
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k\ or mass scale Mi, in order for that region to collapse by 
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[HI- 

Anticipating the linear halo bias that we will calculate 
next, in FIG.|3]we plot the ratio of the flat barrier solu- 
tion, Eq. ([To]), to the full curved barrier solution, Eq. ([9| 
using the numerical procedure of Sec. II F for four halo 
masses. The amplitude of the flat barrier is taken to be 
the collapse threshold Si at high ki, which in the case 
of our models is very close to the Newtonian value. The 
flat barrier approximation is very good for galaxy-mass 
haloes — for a halo of 10^'^ M©, the fractional discrep- 
ancy in the bias is of order 10""^ for the parameters con- 
sidered. Decreasing \a\, increasing A, or increasing Zform 
will flatten the barrier, further reducing this difference. 
Effects due to barrier curvature are more readily appar- 
ent for cluster-mass haloes, because the barriers in FIG. [2] 
are only curved at large masses. 



Formation bias 



a range of parameter values for a and A. At smaller 
ki in the weaker-gravity case (a > 0), the barrier is 
higher because these larger-scale perturbations enter the 
weakened-gravity regime earlier and experience more of a 
reduction in the gravitational force, and so need a larger 
initial overdensity to push them over the edge into col- 
lapse; the converse holds for the stronger-gravity (a < 0) 
model. At high fci the barrier approaches the GR re- 
sult, which can be thought of as the value 1.671, scaled 
back to Zi — 100 using the ACDM growth factor ratio 
D{zi)/D{zfoTm)- While not a focus of this paper, FIG. [2] 
suggests the mass function at the high mass end would be 
significantly affected by the curving of the barrier, espe- 
cially its lowering at high mass for a < 0. The observed 
cluster count already constrains a to be not too negative 



There are two ways in which the halo formation bias 
(i.e. halo bias with Zform = ^obs) differs from the stan- 
dard expression, due to for instance a modification to the 
gravitational law. One is from the scale-dependent lin- 
ear growth, the other is from the curved barrier. The 
former is the dominant effect. As we have seen already 
from FIG. |3j the curvature of the barrier does not have 
a significant impact on the halo bias, on scales fco ^ 0-05 
h/Mpc. The impact on smaller scales is higher, but still 
minor for sufficiently light haloes. 

Our most general expression for the halo bias is given 
by Eq. ([9|. The factor of D{So, z,) / D{So, Zobs) encapsu- 
lates the first effect. For instance, in the Yukawa model 
with a > (weaker gravity on larger scales), this factor 
is larger on large scales, and the halo bias is correspond- 
ingly higher. Physically, what drives the higher halo bias 
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on larger scales is this: recall that we are interested in 
how the halo overdensity 6fi relates to the matter over- 
density S at redshift Zobs; for a fixed S, the overdensity Sq 
at the initial redshift Zi (where the random walk is per- 
formed) is higher on larger scales if a > 0. This larger Sq 
boosts halo formation, and results in a higher halo bias. 

The curved barrier is in principle another way in which 
a scale dependence is imprinted on the halo bias. For 
this to occur, the barrier must absorb some trajectories 
before they reach {Si,6i{Si)) in a way which depends 
on the starting point Sq. We do not observe this in the 
models we consider, because on large scales the random 
walk has a negligible probability of crossing the barrier 
because the variance of each step, S, is so low. Therefore 
large-scale changes to the barrier have very little effect 
on the bias. (If a is strongly negative, say a ~ —0.7, the 
barrier is so low on large scales that it absorbs the random 
walk paths before they reach 'sensible' halo scales, see 
FIG. 



Mo{Ms/h) 
Ifli* 10" 



This produces a strongly scale-dependent bias 
and drastic under-production of low-mass haloes, which 
is clearly incompatible with observations.) 

FIG. [4] shows the halo bias computed using the full 
excursion set theory, with a barrier determined by the 
exact spherical collapse calculation. The scale depen- 
dence of the bias comes almost entirely from the large- 
scale growth factor, and is well described by our flat- 
barrier expressions, Eqs. ([TT|), ([l3|), or ([T5|. As the 



halo mass increases (and so the variance Si decreases) 
the factor multiplying the scale-dependent growth factor 
ratio, [Si'^ / Si — l)/(5i, becomes larger, leading to more 
pronounced scale dependence in the formation bias. The 
different curves in FIG. |4] are produced by the same scale- 
dependent factor multiplying a mass-dependent ampli- 
tude. 
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-0.5, A = lOh-i Mpc. 



C. Passive bias evolution 

In the previous section we observed the haloes at the 
time of their formation. Some galaxy populations form 
at an early redshift and then evolve passively, without 
further major merger events. In this scenario the total 
number of haloes is conserved after the formation epoch; 
the haloes increase their mass only by accreting much 
smaller haloes. The mass assigned to the haloes is their 
mass at formation, disregarding subsequent accretion. 

Here Zform Zobs] the formation redshift determines 
the barrier shape and height via 5i — ^crit,i (this is 
the initial overdensity threshold which has to be evolved 
back from the formation redshift), while the observation 
redshift determines the relationship between 6q and S, 
Eq. ([6]). This combination reproduces, in the excur- 
sion set formalism, the modified-gravity bias-evolution 
results previously derived from the perturbation equa- 
tions by HP [To], in the sense that 5i — 1 scales with 
D{SQ,Z{orm)/D{So,Zohs), as in Eq. ([l2]). However, it 



FIG. 4: Bias versus scale for Zform ~ Zobs = 0. The halo 
masses in each case are, in units of Mq, 5 x 10^'* (most 
biased), lO", 5 x lO", lO", 5 x 10^^ and 10^^ (least biased). 
For the three largest masses, the curves are discontinued be- 
yond the scale at which halo exclusion becomes important, i.e. 
where Mi /Mo = 0.01. The panels illustrate the predictions 
for two Yukawa models. 



general, but is instead given by Eq. (13). 



should be kept in mind that according to the excursion 
set theory, the formation bias is not scale-independent in 



The combined effects of formation bias and passive bias 
evolution can be non-intuitive. For example, in FIG. [5] 
the stronger large-scale gravity (a < 0) causes the bias to 
be smaller on large scales, but passive evolution decreases 
the smaller scale bias faster. This is because on smaller 
scales the bias is greater, and thus is more quickly diluted 
by passive evolution. 
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FIG. 5: Passive bias evolution of a 5 x 10^^ Mq halo 
which forms at Zform = 1 and is observed at Zohs = 1 (top), 
0.5, and (bottom); in this case a = —0.5, A = 10 Mpc. 



D. Comparison with observations — an example 



For comparison with observations, it is sometimes con- 
venient to consider an apparent bias 6app, defined as 



D{ko,Zi) DACBuizohs) 



6i(fco), (20) 



where Dacdm is the growth factor according to the stan- 
dard ACDM model, and D is that according to our non- 
standard model. This is constructed such that the power 
spectrum of haloes in the non-standard model P\ 



halo 



IS 



related to the matter power spectrum in ACDM Pacdm 
by 



haloes 



app 



ik)PA. 



CDM ■ 



(21) 



This is a necessary step when we don't know the mat- 
ter power spectrum directly, and so must compare the 
observed galaxy distribution to what would be expected 
from an underlying standard-gravity matter distribution. 
If the matter power spectrum were known, from weak 
lensing for example, we could compare the measured bias 
to the true bias bi directly. 

The apparent bias is plotted against scale in FIG. |6] 
On large scales the apparent bias has a constant slope, 
which is the same for all halo masses. This is be- 
cause the form of the true bias, 6i w 1 + /i(fco)/2(Mi) 
(Eq. [TTj ) , when multiplied by the scale-dependent factor 
Diko, Zohs )/D{ ko, Zi), becomes &app ~ gi(fco) + 52(^^1)- 

In FIG. 6(c)[ the large-scale apparent bias actually in- 
creases with passive evolution — this is due to the over- 
all multiplicative factor in front of bi in the previous 
equation, which is greater than one on large scales when 
a < 0. The opposite occurs for positive a: the large-scale 
apparent bias drops faster than the physical bias. 

The high intrinsic luminosity of luminous red galax- 
ies (LRGs) makes them useful tracers of the large-scale 
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(c)Passive apparent bias evolution of a 5 X 10^^ Mq halo; 
Zform = l-Oi -^obs = 1-0 (solid) , 0.5 (dashed), (dashed-dotted). 
Compare with FIG. [5] 

FIG. 6: Birth and evolved apparent (rather than physical) 
biases for the same parameters and halo masses as FIG.|4]and 
FIG.H 
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FIG. 7: Galaxy biases derived from the SDSS LRGs, with 
Yukawa and Newtonian best fit curves. The best fit parame- 
ters are a = 0.065, A — 6.5 Mpc 



cosmic density field. We calculate the galaxy bias of the 
LRGs as a function of scale, using the LRG power spec- 
trum derived from the SDSS galaxy redshift survey [iP] . 
and the standard matter power spectrum, whose param- 
eters are set to the WMAP 5-year maximum-likelihood 
values [38]. The LRGs are mostly believed to occur as 
the largest galaxy in individual haloes, and so the galaxy 
bias is approximately the host halo bias. 

A description of the halo clustering with only the first 
bias parameter and the matter power spectrum is con- 
sistent only on linear scales, which we estimate here as 
those scales where the scale dependence introduced into 
the bias by a non-zero iS'o(fco) in Eq. 



10 is small. This 



restricts our analysis of the LRG power spectrum to the 
eleven data points available for /cq ^ 0.06 ft,Mpc~^. 

The LRGs may be a passively evolving population [IT] 
(however, see [12] for an opposing view). Despite this, 
we use the formation bias at the median redshift of the 
sample, Zform = -^obs = 0.35, in order to include any mass 
the haloes may have accreted since they formed. Ide- 
ally, one would bin the galaxies by redshift and account 
for the changing predicted bias shape for each modified 
gravity model, thus providing two tests: the shape of the 
bias curve at each redshift, and the evolution of the bias 
(in both shape and magnitude) with time. This may be 
feasible with future data sets. 

A complete analysis would include marginalization 
over the matter spectral index and tig. To derive some 
preliminary constraints, we fix all the properties of the 
matter power spectrum and (ACDM) expansion history, 
and marginalize over halo mass in a minimization pro- 
cedure^ for each combination of a and A. Fitting the 
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Here, we do not distinguish between the fc's in our equations 
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FIG. 8: Contours of constant Ax^ in the q-A plane. The 
regions containing a = are 'included' by the contours. 



apparent formation bias at the median redshift of the 
LRG sample gives best-fit parameters a — 0.065, A = 
6.5 h"^ Mpc, Ml = 1.2 X lO^^ /j-i Mq, with = 10.08. 
This is statistically indistinguishable from the Newtonian 
fit, with = 10.21. In fact, one could argue that the 
data does not favor the Yukawa model since it requires 
more parameters without improving much. However, 
the error bars are at the moment rather large. The ob- 



(which are really defined in real space, in the sense of fc = 1/_R) 
and the true Fourier-space fc's as in the observed galaxy power 
spectrum. This is not an issue if one has a scale-independent 
growth factor, and therefore a scale-independent large-scale bias. 
But for our present purpose, some extra care should in principle 
be taken in determining the correct mapping between our fc's 
and the observed Fourier-based fc's. Simulations can be used to 
address this issue, which we leave for future work. 



11 



served biases and model fits are shown in FIG. [71 

FIG. [s] shows contours of constant Ax^ in the a~X 
plane. Marginalization over a or A is poorly constrain- 
ing, since in both cases the Newtonian limit {a ^ 
or A — >■ oo) is allowed, and so the marginalized Ax^ 
cannot be larger than 0.13. This method places bet- 
ter limits on stronger-than-Newtonian-gravity models. 
These constraints should be considerably improved by 
future galaxy clustering measurements, especially those 
on scales fco < 0.1/iMpc~^. 



IV. CONCLUSIONS 

We have worked out the large-scale halo bias bi in a 
theory that has a scale-dependent linear growth, such 
as modified gravity or clustered dark energy. We give 
several expressions in |II C| & |II D[ The most general ex- 
pression is Eq. (|9|. It relates bi to the first crossing 
distribution /, which can be computed using the method 
described in §11 F[ for a general curved barrier. Note that 
all quantities describing the random walk: Si, Si,6o, Sq 
are defined at the initial redshift Zi, as opposed to the 
formation redshift Zform- 

A much simpler expression, Eq. ( |15[ ), obtains in cases 
where the barrier can be approximated as flat. As can 
be seen in the example depicted in FIG. [S] the barrier is 
generally not flat once a scale-dependent growth factor 
is allowed. However, it is often true that the barrier 
deviates from flatness only for the largest masses, and a 
fiat barrier approximation (set at the small mass level) 
actually works quite well in predicting the halo bias for 
more modest halo masses (FIG.|3]). The expression in Eq. 
(151 is particularly easy to use because the threshold for 
collapse (5i.form and the variance S'ljorm are defined at 
the formation redshift as usual. 

Our main result is the imprinting of a scale dependence 
on the formation bias, by the modified linear growth fac- 
tor. This can be understood to occur because regions 
of the same overdensity but different sizes today would 
not all have had the same overdensity at an earlier red- 
shift, and hence would have varying halo densities (see 
Sec. IIIB). The scale dependence of the growth factor 



manifests itself at a single redshift, making it possible 
to look for it with a local-universe snapshot. Using this 
effect and an appropriate galaxy sample, one can place 
constraints on gravity theories, as we demonstrate with 
a simple Yukawa theory. 

We have focused on the case in which the formation 
bias is deterministic. However, this bias may be stochas- 
tic, in which case this stochasticity will also be scale- 
dependent [1P|. Although we will not perform a detailed 
analysis of stochastic bias here, we note that, even if the 
formation bias is deterministic in, say, fc-space, the scale 
dependence means that it will be stochastic in real space 
(and vice- versa) [l^. Thus, the stochasticity of the bias 
is another potential signature of modified gravity theo- 
ries. 



Looking forward, it would be interesting to work out 
how a scale-dependent large-scale bias can be incorpo- 
rated into optimal weighting schemes for recovering the 
mass distribution from a galaxy or halo catalog [JH US] . 
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Appendix A: Density field extrapolation 

In an excursion set calculation, there are in principle 
four different redshifts that might be relevant: the initial 
redshift Zi (taken to be 100 in this paper), the forma- 
tion redshift Zform, the observation redshift Zobs, and the 
redshift at which all the random walk variables (the vari- 
ances So, Si, etc) are defined: let us call it ^ficidsj or the 
random walk redshift. 

In the standard theory, i.e. GR, the choice of Zficids 
is immaterial. The simplest way to see this is to exam- 
ine the flat barrier (but otherwise general) expression for 
bias in Eqs. (12) and (13 1. In this expression, the ran- 



dom walk variables Si and 5*1 are defined at the initial 
redshift Zi. In other words, the choice Zfioids = Zi has 
been made. It is simple to see that if the growth factor 
D were scale-independent, choosing a different random 
walk redshift, i.e. Zi — >■ Zfloids 7^ Zi would have made no 
difference to the predicted bias. This is because in the 
expression for the formation bias (Eq. 



13 ), Sl/Si 



dependent of the random walk redshift by definition, and 
6iD{So,Z{orm)/D{So ) (with Si defined at Zfidds) is 

also independent of Zfjdds- The key is that the growth 
factor relevant for 5i, which is on scale Si, is the same 
as the growth factor on scale 5*0. This is true in GR, but 
not true for non-standard theories where the growth fac- 
tor is generally scale-dependent. Indeed, in the standard 
theory, it is common practice to choose Zfieids = ^form, 
such that for instance in an Einstein-de Sitter universe, 
5i takes the value 1.686. 

In non-standard theories, the formation bias given 



by Eq. ( 13 ) is still insensitive to the precise choice of 
•Zficids — Zi, as long as Zi is sufficiently early. This is 
because we expect the (sub-Hubble) growth to be scale- 
independent at sufficiently early times. Choosing Zfieids 
to be some late redshift, however, would change the pre- 
dicted bias. Our choice of Zfidds being equal to some 
early initial redshift is in keeping with the excursion 
set philosophy, namely that characteristics of the initial 
random Gaussian fiuctuations, such as how they change 
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with smoothing scale (i.e. the random walk), determine 
whether a region ultimately collapses into a halo by a 
certain rcdshift. 



Appendix B: Scalar-tensor Theory 

In this section, we would like to make explicit how the 
Yukawa model we use as an illustrative example can be 
obtained from a scalar-tensor theory . In the Einstein 
frame, a scalar-tensor theory with potential scalar self- 
interactions can be written as 



1 



SttG 



1 



:R 



V,,</7V^<^ - V{ip) 



(Bl) 



where R is the Ricci scalar, the dimensionless (p medi- 
ates a scalar fifth force, V{ip) is the scalar potential, and 
1pm represents matter which is minimally coupled not to 
the Einstein frame metric g^jy but its conformal cousin 
ri~^((/?)gpj,. Here, ip is typically small such that the con- 
formal factor can be expanded as 



n^{(p) - 1 - 2/?(^ 



(B2) 



where /3 is a coupling constant, typically of order unity 
for a gravitational-strength scalar force. 

The total gravitational -I- scalar acceleration an in- 
finitesimal particle experiences is 



X, 



-dicj) = ~di(j) - I3di(p , 



(B3) 



where </) is the effective (total) gravitational potential, 
and (f) is the (Einstein frame) gravitational potential. The 
scalar field tp and the gravitational potential (j) satisfy the 
equations: 



dV 

dp 



dV_ 

d(p 



Sn^Gpa^S (B4) 



where S is the overdensity, p is the mean density, and (p is 
the scalar field value at mean density. If the potential V 
were negligible, one can see that (f) and ip would be pro- 
portional to each other, and the test particle would accel- 
erate according to an effective gravitational potential of 
(j) = (j>{l + 2l3'^). More generally, assuming the fluctuation 
of ip from (p is small, one can Taylor expand and rewrite 
the potential terms on the right as ~ {a^ /\^){p — ip), 
where 1/A^ is the second derivative of the potential eval- 
uated at (p. One therefore finds ip[l+a? / (fc^A^)]/(2/3) = 
in Fourier space. This means the total effective gravita- 
tional potential </> = [l-f 2l3^P/{k'^ + a^/A^)] . Com- 
bining this with the Poisson equation for </> gives us the 
Yukawa model in Eq. (18), if one makes the identifica- 
tion a = 2/32, and Gat = G(l 2/3^). The Yukawa model 
we study assumes A is a constant, but it should evolve in 



general, depending on the precise form of the potential 
V. 

Note that a is always positive in this sort of formu- 
lation. The case of a negative a is best regarded as 
purely phenomenological. Note also that depending on 
the potential V, the scalar-tensor theory could exhibit a 
chameleon mechanism on sufficiently small scales, under 
which [ip — ip)/ip cannot be considered a small perturba- 
tion. In this case, one finds 3 regimes: the largest scales 
(compared to A) on which the scalar force is Yukawa sup- 
pressed, the intermediate scales on which the scalar force 
is operative, and the small scales on which the chameleon 
mechanism takes over to suppress the scalar force. With 
the chameleon mechanism, there is the additional pos- 
sibility of equivalence principle violations, namely that 
macroscopic objects (as opposed to test particles) do not 
necessarily all fall at the same rate. See H?! US] for fur- 
ther discussions (note that their a is the same as (3 here) . 



Appendix C: Spherical collapse 



The Euler equation in an expanding universe is 



V)v a' „ 

— = v-Vq 

Drj a 



(CI) 



where 77 is conformal time, v is the peculiar velocity 
da;/d77 [x being comoving position), a' = da/drj, V is 
the gradient operator with respect to the comoving co- 
ordinates, and 4> is the perturbation in the gravitational 
potential. D/D77 is strictly speaking a material deriva- 
tive, but, since we are going to follow the motion of a 
surface (a spherical shell), we are using Lagrangian co- 
ordinates, and so can think of this as just a normal time 
derivative d/dry. Rewriting this in terms of the physical 
(proper) co-ordinate i? of a spherical shell, with the scale 
factor as our time co-ordinate, we find 



da2 



1 



H 



dH 
aa 



R-a 



dR 
da 



(C2) 

We investigate three approaches to calculate the evo- 
lution of a spherically-symmetric perturbation, with an 
initial top-hat profile. In the first case, we divide the per- 
turbation into spherical shells and follow the motion of 
each shell under the gravitational force from every other 
shell. This method includes the effect of inhomogene- 
ity (i.e. departure from a pure top-hat profile), which 
we might expect to be important because of the viola- 
tion of Gauss's law. The second method approximates 
the perturbation as a homogeneous sphere at every time 
step, and only solves for the motion of the outermost 
surface. T hird ly, we use Raychaudhuri's equation to re- 
place Eq. ( C2 1 , approximating the problem as local. In 



all cases we use a fourth-order Runge-Kutta integrator 
with Cash-Karp co-efficients and adaptive step-size ad- 
justment, similar to that found in Press et al. [35] . 



13 



The initial conditions for the velocity are 



dR 
da 



R 
a 



1 



3(1 + 5) 



dS 
da 



the same as would be used in Newtonian gravity, because 
these are set early in the matter era, when the effects of 
any modifications to gravity will be negligible. 



a. Inhomogeneous simulations 

For the inhomogeneous simulations, the force experi- 
enced by any shell can be split into three contributions: 



that due to the Newtonian l/r part of Eq. (19 1, and 



those due to the additional Yukawa-type contributions 
from shells internal and external to the point of inter- 
est. Integrating the potential measured at r sourced by a 
shell of thickness dr' centered on r', and taking the gra- 
dient with respect to r, we find that the Newtonian-type 
contribution is 



Gn M,nt-{4Tr/3)p{ay 
I + a 



(C3) 



where Mi^t is the total mass internal to the proper radius 
r; the Yukawa-type contribution from internal shells is 

^ 1 ^ Gn , r'\ (\ 1 

Vr0Y,int = 2™— — dr p[a)6i - + t 

\ -V a r \r A 

X {exp [-{r - r')/A] - exp [-(r + r')/A]} , 

(C4) 

where Si is the overdensity of the source shell, and the 
Yukawa-type contribution from external shells is 

, G ]\[ I _, , r'A 

Vr'/)Y,cxt = -2-Ka— — dr /9(a)dj — 
1 -I- a r 



exp [— (r' — r)/A] 




+ - + - cxp[-(r' + r)/A] 



(C5) 



h. Homogeneous simulations 

For the homogeneous-sphere simulations, the gradient 
of the gravitational potential at the surface of a pertur- 
bation of total mass M and radius R is 



G 



N 



M ■ 



47r 



p{a)R: 



1 

i?2 



R 



1 + a ( 1 + - ) F{R/X) 



(C6) 



and the form factor [49] is given by 
3 



rC ^ {a;cosh(a;) — sinh(2;)} 




0.6 0.8 l.( 
Scale factor, a 



FIG. 9: Comparison of nonlinear solution methods. Lower 
solid line: Newtonian solution; upper solid line: homogeneous 
simulation; dashed line: inhomogeneous simulation; dotted 
line: local transformation (evaluated at center). 



c. Local density-field transformation 

Gaztanaga & Lobo [SD] use a combination of the conti- 
nuity and Raychaudhuri equations to derive a differential 
equation for the nonlinear evolution of the density field 
in Brans-Dicke gravity. We adopt their method to our 
problem. 



Taking the divergence of the Euler equation (CI I, and 



ignoring shear and vorticity, it can be shown that 



DO 



1 



(C7) 



where it should be kept in mind that D/Drj is a La- 
grangian derivative, and 9 = 'V - v. Putting this together 
with the continuity equation DS/Dr] = —(1 + 6)9, we 
obtain: 



a' DS 
a Drj 



4 1 DS 
?,l + S~Dq 



= (1 + 5)V2(/. (C8) 



Note that spatial derivatives in this equation are comov- 
ing. To close this equation, we need to relate V^</) to 



S. Taking the divergence of Eq. (19 1, and applying the 
divergence theorem, we find located at the center 
of the spherical top hat (where shear and vorticity can 
justifiably be ignored) to be: 



. Gn _ 2r 

47r pa 

1 + a 



1 



1 



-R/\ 



(C9) 



The proper radius R and S are related by mass conser- 
vation, i.e. (1 -I- S)R^ /a^ is a constant. 



d. Comparison of methods 

The three modified-dynamics methods give very sim- 
ilar results when the deviation from the Newtonian tra- 
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jectory is small, usually when i? < A at all times. How- 
ever, in more extreme scenarios the Raychaudhuri equa- 
tion method overstates the degree to which the gravity 
modification changes the trajectory, in comparison to the 
inhomogeneous sphere method, which we assume is the 
most accurate. R{a) trajectories for the three modified 
gravity solutions are shown in FIG. [9j along with the 



associated Newtonian path; a = 5, X = 5 h ^ Mpc. 

It appears that the local approximation inherent in 
the Raychaudhuri method breaks down when R > X. 
However, the homogeneous-sphere approximation con- 
tinues to hold until well beyond the limits of interest: 
M < IQis Mg/h, A > 10 h-^ Mpc. 
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